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Abstract We overview nonequilibrium Green function combined with density functional 
theory (NEGF-DFT) modeling of independent electron and phonon transport in nanojunc- 
tions with applications focused on a new class of thermoelectric devices where a single 
molecule is attached to two metallic zigzag graphene nanoribbons (ZGNRs) via highly trans- 
parent contacts. Such contacts make possible injection of evanescent wavefunctions from 
ZGNRs, so that their overlap within the molecular region generates a peak in the electronic 



^ . transmission. Additionally, the spatial symmetry properties of the transverse propagating 

^ ' states in the ZGNR electrodes suppress hole-like contributions to the thermopower. Thus op- 

timized thermopower, together with diminished phonon conductance through a ZGNR | molecule | ZGNR 
inhomogeneous structure, yields the thermoelectric figure of merit ZT ~ 0.5 at room tem- 
perature and 0.5 < ZT < 2.5 below liquid nitrogen temperature. The reliance on evanescent 
mode transport and symmetry of propagating states in the electrodes makes the electronic- 
transport-determined power factor in this class of devices largely insensitive to the type 
. of sufficiently short conjugated organic molecule, which we demonstrate by showing that 

Q ' both 18-annulene and CIO molecule sandwiched by the two ZGNR electrodes yield sim- 

O . ilar thermopower. Thus, one can search for molecules that will further reduce the phonon 

thermal conductance (in the denominator of ZT) while keeping the electronic power fac- 
tor (in the nominator of ZT) optimized. We also show how often employed Brenner em- 
^ ' pirical interatomic potential for hydrocarbon systems fails to describe phonon transport in 

I our single-molecule nanojunctions when contrasted with first-principles results obtained via 

O '• NEGF-DFT methodology. 
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1 Introduction 

1.1 Why study nanoscale thermoelectrics? 

Thermoelectrics transform temperature gradients into electric voltage and vice versa. Al- 
though a plethora of widespread applications has been envisioned, their usage is presently 
limited by their small efficiency (I]. Thus, careful tradeoffs are required to optimize the 
dimensionless figure of merit 

ZT = S^GT/k, (1) 

which quantifies the maximum efficiency of a thermoelectric cycle conversion in the linear- 
response regime. This is due to the fact that ZT contains unfavorable combination of the 
thermopower S, average temperature T, electrical conductance G and thermal conductance 
K = Kf.\ + ifph. The total thermal conductance has contributions from both electrons (Qi and 
phonons (Cpj,. The devices with ZT > 1 are regarded as good thermoelectrics, but values of 
ZT > 3 are required for thermoelectric devices to compete in efficiency with mechanical 
power generation and refrigeration ITl . 

The traditional efforts to increase ZT have been directed toward selective reduction of 
the lattice thermal conductivity (Cph, using either complex (through disorder in the unit cell) 
bulk materials O or bulk nanostructured materials (Sj, while at the same time maintaining 
as optimal as possible electronic properties encoded in the power factor S^G. However, 
decennia of intense research along these lines have increased ZT of bulk materials only 
marginally |T||. A complementary approach engineers electronic density of states to obtain a 
sharp singularity (3]|4|| near the Fermi energy Ep which can enhance the power factor S^G, 
such as in Tl-doped PbTe where ZT ~ 1.5 has been reached at 775 K 0. 

The nanoscale and low-dimensional (6) devices offer additional degrees of freedom that 
can be tailored to achieve high-Zr, as exemplified by the recent ground-breaking exper- 
iments demonstrating how rough silicon nanowires (SiNW) can act as efficient thermo- 
electrics {ZT ~ 0.6 at r = 300 K) although bulk silicon (ZT = 0.01 at 7 = 300 K) is not (7] 
[8). Another example of nanoscale thermoelectrics has emerged from the recent experiments 
measuring thermopower of single-molecule nanojunctions t9l fT0lllllll2lll3ll and quantum 
dots fT4ll . 

Availability of efficient nanoscale thermoelectrics could make possible targeted cooling 
of local hotspots I15II due to the ease of on-chip integration. To make use of low-dimensional 
thermoelectric devices in macroscale applications, many nanowires must be placed in par- 
allel, so issues of the nanowire size and packing density arise (6). 

Besides device applications, the search for optimal ZT has ignited basic research in 
condensed matter physics and various engineering disciplines aimed at deepening our un- 
derstanding of heat flow in nanoscale or in unconventional bulk systems. For example, the 
recent review article (161 on heat flow and thermoelectricity in single-molecule nanojunc- 
tions and atomic wires highlights that even apparently basic issues are not well understood 
in such systems. In conventional systems, where heat and charge currents are transported 
by Landau quasparticles, ZT is normally limited by the Wiedemann-Franz law stating that 
Ke\/GT is a system-independent constantQ However, the Wiedemann-Franz law is a result 
of Fermi-liquid theory and breaks down 1181 in correlated bulk materials or in nanoscale 

' We should mention here that the Lorenz ratio K^i/GT calculated for transport of non-interacting 
electrons through several single-molecule nanojunctions shows variations by tens of percent from the 
Wiedemann-Franz law as the chemical potential crosses a transmission resonance, and much larger devia- 
tion around the transmission nodes 1171 . 
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systems (such as quantum dots or metallic islands) with strong Coulomb interaction effects. 
This has necessitated the development of novel theoretical techniques 1191 to tackle thermo- 
electricity in correlated bulk materials, such as the combination of the dynamical mean field 
theory 1201 with the local density approximation and the Kubo formula, which has revealed 
enhanced thermopower in, e.g., FeSb2 and NatCo02 due to electronic correlations |21I . In 
the realm of nanoscale correlated systems, significant ZT values were predicted for, e.g., 
Kondo con^elated quantum dots I22L metallic single-electron transistors 1181 and Kondo 
insulator nanowires 11231 . 



1.2 What is interesting about thermoelectricity in single-molecule nanojunctions? 

Very recent experiments l9l ll0lfTT]|12lll3l have measured thermopower S as induced thermo- 
electric voltage in response to a temperature difference across organic molecule sandwiched 
between two gold electrodes. This has ignited vigorous theoretical efforts 111611 17ll2"4ll25ll26i 
I27]|28ll29|[30ll31ll32ll33ll34ll35ll36ll37ll38l to explore devices where a single organic molecule 
is attached to metallic or semiconducting I31II electrodes. In such nanojunctions, the di- 
mensionality reduction and possibly strong electronic correlations |'34|35|361 allow for the 
increase in S concurrently with diminishing (Cph while keeping the nanodevice disorder- 
free 1391 . For example, creation of sharp transmission resonances near the Fermi energy 
Ep by tuning the chemical properties of the molecule and molecule-electrode contact via 
chemical functionalization can substantially enhance II26II31I the thermopower S which de- 
pends on the derivative of the conductance near Ep. At the same time, the mismatch in the 
phonon density of states at the electrode | molecule interface can severely disrupt phonon 
transport 1401 . when compared to homogenous clean electrode made of the same material, 
thereby leading to small (Cph. 

Besides offering a new route toward high-ZT devices, thermoelectric properties of single- 
molecule nanojunctions have been investigated as a tool that could resolve a number of fun- 
damental issues in molecular electronics |41S. For example, Ref. 1421 has suggested that 
thermoelectric voltage in single-molecule nanojunctions, which is large enough and rather 
insensitive to the details of coupling to the contacts, can be employed to locate the position of 
the Fermi energy relative to the molecular levels. This concept has sparked the development 
of experimental techniques that can unambiguously identify the molecular orbital closest 
to the Fermi level of the electrodes by measuring thermopower and current-voltage char- 
acteristics of single-molecule nanojunctions 1131 . The measurement of S could also probe 
room-temperature quantum interference effects in transport through single-molecule nano- 
junctions I34II43 Il44ll45 Il46p7ll where Ref. (175 has predicted dramatic enhancement of the 
thermopower near a transmission node because the flow of entropy (as an inherently inco- 
herent quantity) is not blocked by destructive quantum interference which block the charge 
current and lead to vanishing conductance. 

The single-molecule nanojunctions also provide access to nonlinear thermoelectric prop- 
erties. Most theoretical studies of nanoscale thermoelectrics have focused on the linear- 
response regime (i.e., close to equilibrium) where one operates close to the small voltage 
V = —SAT which exactly cancels the current induced by the small temperature bias AT . 
As ZT — )• o°, the efficiency approaches the ideal Camot value 7] — )■ 7]c = 1 ~ T /{T + AT). 
However, in the linear-response regime AT <^T typically investigated for bulk materials, 
the efficiency stays low rj k AT /T even if ZT can be made very large. The recent experi- 
ments |9l on single-molecule nanojunctions showed nonlinearities in the measured S already 
at zir w O.ir. Furthermore, the analysis of simple phenomenological models of nanojunc- 
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Fig. 1 Schematic view of (a) ZGNR| 18-annulene|ZGNR and (b) ZGNR|C10|ZGNR single-molecule nano- 
junctions. The contact between the source and the drain 8-ZGNR (consisting of eight zigzag chains) metallic 
electrodes and molecules is made via five-membered rings of carbon atoms (dark blue). The hydrogen atoms 
(light yellow) are included to passivate the edge carbon atoms in the nanoribbon or danghng bonds in the 
molecules. 

tions has suggested I48II36I that optimal thermoelectric operation can be achieved in the 
out-of-equilibrium nonlinear regime. Such regime also requires theoretical approaches that 
go beyond usual concepts of the figure of merit ZT and thermopower S defined in the linear- 
response regime 1361 . 



1.3 What is the scope of this mini-review? 

This article has two principal goals: (/) to present some technical details on the usage 
of nonequilibrium Green function combined with density functional theory (NEGF-DFT) 
computational methodology to study electron and phonon transport quantities which are re- 
quired to understand thermoelectricity in single-molecule nanoj unctions; (n) to overview 
applications of NEGF-DFT methodology to a new class of single-molecule-based thermo- 
electrics 1321 illustrated in Fig.[T] 

In these devices, organic molecule is attached to metallic graphene nanoribbons with 
zigzag edges (ZGNR) via strong covalent bond which provides high transparency of ZGNR| molecule 
contact. This allows for evanescent wavefunctions to penetrate from the electrodes into the 
molecular region thereby generating a transmission resonance which, together with spatial 
symmetry of the transverse propagating states within the ZGNR electrodes, yields highly 
optimized power factor S^G. Using two different molecules in Fig.[T] we demonstrate that 
S^G optimized in this fashion is independent on the type of the molecule employed, which is 
a feature that can be exploited to search for the molecule that will bring the largest reduction 
of the phonon thermal conductance. 
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The article is organized as follows. In Sec.|2]we discuss NEGF-DFT approach to the 
computation of electronic conductance, thermopower and thermal conductance for realistic 
single-molecule nanojunctions built from carbon and hydrogen atoms. The same approach 
applied to phonon thermal conductance is discussed in Sec.[3]for the elastic regime where 
electron-phonon I49II50I interactions or phonon-phonon 1511 scattering processes are ne- 
glected. The electronic transmission function and thermopower for ZGNR| 18-annulene|ZGNR 
and ZGNR|C10|ZGNR junctions are discussed and compared in Sec.|4] In Sec.|5]we show 
how the usage of ZGNR electrodes, whose transverse propagating eigenstates have spe- 
cial symmetry properties in real space, makes it possible to block the hole-like portion of 
the transmission function thereby enhancing the thermopower. Section[6]shows the phonon 
transmission function and the corresponding thermal conductance for these two junctions. In 
Sec.[6l we also compare phonon transport quantities obtained via full NEGF-DFT method- 
ology to computationally faster NEGF coupled to Brenner empirical interatomic potential 
(EIP) where we find failure of the latter technique to describe nanojunctions in Fig. [T] Us- 
ing the quantities computed in Sec. |4] and Sec.|6l we construct the thermoelectric figure of 
merit ZT in Eq. (|T]| for the two single-molecule nanojunctions in Sec.|7] We summarize and 
discuss briefly future directions for the computational modeling of single-molecule thermo- 
electrics in Sec.[8l 

2 NEGF-DFT methodology for the computation of electron transport-determined 
thermoelectric quantities 

The traditional tools for the computation of thermoelectric transport coefficients, such as 
the semiclassical Boltzmann equation 1521 . cannot be used for quasiballistic nanometer- 
size active region attached to much larger reservoirs The proper description of such open 
quantum system can be achieved using quantum master equations for the reduced density 
matrix of the active region I54ll55ll56ll57l or NEGF formalism I58II . The former is typically 
used when the active region is weakly coupled to the reservoirs (so that coupling between 
the molecule and the electrodes is treated perturbatively), while the latter is employed in 
the opposite limit I59II . The complementary nature of the two methods and the boundaries 
of their validity were analyzed in, e.g., Ref. 1601 using simple phenomenological models of 
single-molecule nanojunctions. 

The NEGF formalism for steady-state transport operates with two central quantities, the 
retarded G(£') and the lesser Green functions (£), which describe the density of available 
quantum states and how electrons occupy those states, respectively I58II . Its application to 
electronic transport is often combined 1611 with the tight-binding (TB) Hamiltonian whose 
hopping parameters are fitted using more microscopic theory!^ This procedure is exemplified 
by the recent calculations 1641 predicting ZT ~ 3 for sufficiently long SiNW of 2 nm diam- 
eter with surface disorder, where electronic subsystem in silicon was described by sp-'d^s* 
TB Hamiltonian with nearest-neighbor hoppings and ten orbitals per site. 

However, such usage of NEGF formalism, where the device Hamiltonian is known 
from the outset, is not suitable for the description of realistic single-molecule nanojunc- 
tions where organic molecule consists of carbon and other atomic species, or in the case 
of GNR electrodes where different atoms or atomic groups are used to passivate dangling 

^ For comparison between Boltzmann semiclassical and Landauer quantum transport approaches applied 
to thermoelectric transport coefficients of conventional translationally invariant systems see Ref. 1531 . 

^ For example, in the case of either bulk graphene 1621 or GNRs 1631 one has to employ TB Hamiltonian 
with up to third-nearest-neighbor hoppings in order to match the DFT-computed band structure. 
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bonds along the edges I65II . In such cases, first-principles input about atomistic and elec- 
tronic structure is necessary in order to capture charge transfer between different atoms 
in equilibrium, geometrically-optimized atomic positions of the molecular bridge including 
molecule-electrode separation in equilibrium, and forces on atoms when they are perturbed 
out of equilibrium. For example, Ref. [6611 shows how linear-response conductance of GNR- 
based devices is computed incorrectly if charge transfer between edge hydrogen and interior 
carbon atoms is not taken into account. 

The state-of-the-art approach that can capture these effects, as long as the coupling be- 
tween the molecule and the electrodes is strong enough to ensure transparent contact and 
diminish Coulomb blockade effects II35II67I . is NEGF-DFT. The DFT part of this frame- 
work is employed using typical approximations (such as LDA, GGA, or B3LYP) for its 
exchange-correlation functional (68|. The sophisticated computational algorithms 16611671 
169)170117111721 developed to implement the NEGF-DFT framework over the past decade can 
be encapsulated by the iterative self-consistent loop: 

;7'"(r) => DFT^ HKs[«(r)] => NEGF ^ n°"'(r). (2) 

The loop starts from the initial input electron density n'"(r) and then employs some standard 
DFT code 11681 typically in the basis set of finite-range orbitals for the valence electrons 
which allows for faster numerics and unambiguous partitioning of the system into the active 
region and the semi-infinite ideal electrodes. The DFT part of the calculation yields the 
single particle Kohn-Sham Hamiltonian 

^Ks[«(r)] = -^+y^«(r), (3) 

y^ff(r) = y//(r) +Ke(r) +yext(r). (4) 

Here ^^^(r) is the DFT mean-field potential due to other electrons where Vnir) is the 
Hartree, Vxclr) is the exchange-correlation and Vext(r) is the external potential contribution. 
The inversion of 7fKs[w(r)] yields the retarded Green function G(£') [see Eq. ((S} below] 
whose integration over energy determines the density matrix via 

P=^J dEG<{E) = J dEG{E)[ifL{E)rL{E) +ifR{E)rR{E)]G\E). (5) 

Here the coherent transport regime (i.e., electron-phonon or electron-electron dephasing 
processes are absent) is assumed, so that G^{E) can be expressed in terms of G(£'). The 
matrix elements ;7°"'(r) = (r|p|r) are the new electron density as the starting point of the 
next iteration. This procedure is repeated until the convergence criterion 1 1 — p™\ \ < 5 
is reached, where 5 <C 1 is a tolerance parameter. The efficient computational algorithms 
for the construction of the density matrix in Eq. ([S} for two- terminal device are discussed in 
Refs. I66II70L and the recently developed algorithms for p in multiterminal devices (includ- 
ing multiterminal thermoelectric s 021 ) can be found in Ref. 11731 . 

The representation of the retarded Green function in the local orbiteil basis {(^, } requires 
to compute the inverse matrix 

G{E) = [ES-n-EL(,E)-ER{E)]-\ (6) 

where the Hamiltonian matrix H is composed of elements Hij = {(j)i\H^s\'l>j)- The overlap 
matrix S has elements Sjj = The non-Hermitian matrices Si^i{[E) represent the re- 

tarded self-energies due to the "interaction" with the left (L) and the right (R) electrodes 1581 
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ICTI . These self-energies determine escape rates of electrons from the active region into the 
semi-infinite ideal electrodes, so that an open quantum system can be viewed as being de- 
scribed by the non-Hermitian Hamiltonian Hopen = Hks [n(r)] -|- ^l{E) + Sr{E). The ma- 
trices riii{E) = i[Si^ii{E) — Slif{E)\ account for the level broadening due to the coupling 
to the electrodes I58II61I . 

The retarded Green function is required to post-process the result of the DFT loop and 
obtain the transport quantities. For example, the zero-bias electron transmission function 
in the elastic (no electron-electron I34II or electron-phonon I49II50I interactions) transport 
regime between the left and the right electrode is given by: 

^,i{E) =Tr{rR{E)GiE)rL{E)G\E)} . (7) 

Each electrode terminates at infinity into macroscopic reservoirs where electrons are as- 
sumed to be thermalized with their Fermi distribution function being f{E — ^i) = fiiE) or 
f{E — /Ir) = fR{E) for the left and right reservoirs, respectively. The transmission function 
■%i{E) in Eq. ^ is then employed to compute the following auxiliary integrals I74II 

K„{^i) = ^l dE%m{E-fir(-^^^), (8) 



where f{E,^) = {l +exp[(Z: — ^)/kBT]}^^ is the Fermi distribution function at the chem- 
ical potential jU (which is the same for both reservoirs in the linear-response regime fii — 
0). The knowledge of K„{fi) finally yields all electronic quantities 



e'-Koin), (9) 

(10) 



^ ^2(M)-[^i(M)]Vgo(M) _ (^j) 

that enter into the expression for ZT in Eq. (|T}. 

Our MT-NEGF-DFT code I451I73I , employed to obtain results in Sec.|4]via the formulas 
discussed in this Section, utilizes ultrasoft pseudopotentials and Perdew-Burke-Ernzerhof 
(PBE) parametrization of the generalized gradient approximation (GGA) for the exchange- 
correlation functional of DFT. The localized basis set for DFT calculations is constructed 
from atom-centered orbitals (six per C atom and four per H atom with atomic radius 8.0 
Bohr) that are optimized variationally for the electrodes and the central molecule separately 
while their electronic structure is obtained concurrently. 

We note that the well-known tendency of DFT-PBE to underestimate HOMO-LUMO 
((HOMO-highest occupied molecular orbital; LUMO-lowest unoccupied molecular orbital) 
energy gap in molecules, and hence overestimate the conductance ||75I does not influence 
the results presented in Sec.|4] This is because for the nanojunctions studied here, the trans- 
mission around the Fermi level is not determined by the off-resonant tunneling through 
the HOMO-LUMO gap. Rather electrons propagate via resonant evanescent modes of the 
graphene electrodes. The energy and shape of such wavefunctions are determined by the 
electronic structure of the graphene nanoribbons and should be well described by DFT. 
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3 NEGF-DFT methodology for the computation of phonon transport-determined 
thermoelectric quantities 

In contrast to rapid development I58II61II of NEGF-TB and NEGF-DFT-based quantum 
transport methods for electron propagation through nanostructures discussed in Sec.[2l com- 
parable methods for phonon transport have emerged relatively slow 1761 . For example, clas- 
sical molecular dynamics (MD) and the Boltzmann equation I77II are widely used traditional 
methods in phonon transport. However, the MD method I78II is not accurate below the De- 
bye temperature, while the Boltzmann equation cannot be used in nanostructures without 
translational invariance. In both cases, the quantum effects become important. 

Only recently, the NEGF formalism has been extended to study the quantum phononic 
transport 1761 . Nevertheless, among recent theoretical studies of single-molecule thermo- 
electric devices using first-principles quantum transport frameworks II24II25||261|27II28II29I 
[30ll31ll32l[33l , most have focused on computing the thermopower S, with only a few 1281 
|3T1|32||331 employing DFT to obtain forces on displaced atoms and then compute (Cph from 
first-principles. The experimental data on the thermal conductance of single-molecule nano- 
junctions is even more scarce I12II79I . 

The phonon thermal conductance in the absence of phonon-phonon 1511 or electron- 
phonon 1149 II 5 01 scattering is obtained from the phonon transmission function J^pii{co) using 
the Landauer-type formula 1801 : 

fj2 ~ ^ ghco/kBT 

'^P" = ^coa>-^pH(c.) _ (12) 



The phonon transmission function ^h(ft>) in such elastic transport regime can be expressed 
in terms of NEGFs for the active region (molecule -I- portion of electrodes) attached to two 
semi-infinite electrodes 

,5'ph(a)) = Tr{ylR(a))D(a))Az.(Z?)D^(ffl)}, (13) 

where the phonon GF is obtained in the same fashion as the electronic one in Eq. l|6ll but 
with substitutions H K, £8 -> w~M and X'l r i7£_«: 

D= [ffl^M-K-TTi-iTfi]"'. (14) 

Here K is the force constant matrix, M is a diagonal matrix with the atomic masses, Ui r 
are the self-energies, and Air{co) = i[Ai f;{co) — A]^^(a))]. 

We compute the force constant matrix K using GPAW (SD, which is a real space elec- 
tronic structure code based on the projector augmented wave method 1821 . The electronic 
wavefunctions are expanded in atomic orbitals with a single-zeta polarized basis set, and 
PBE exchange-correlation functional is used. The whole active region, which includes 27 
layers of ZGNR electrodes, is first relaxed to a maximum force of 0.01 eV/A per atom. 

Subsequently, we displace each atom I by Q/a in the direction a = {x,y,z} to get the 
forces Fjjj {Qia) on atom 7 7 in direction jS. The elements of K-matrix are then computed 
from finite differences 

Fjl}{Qla)-Fjp{-Q,a) 
J^Ia.Jp = ^71 ■ (Ij) 

The intra-atomic elements are calculated by imposing momentum conservation, such that 
Kiajp = —^j^iKjajp- 
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The methodology described in this Section does not take into account resistive umklapp 
phonon-phonon scattering which plays an important role in interpretation of experiments on 
room temperature thermal conductivity of large-area graphene ||83H84l l85|. However, this ef- 
fect, which is easy to describe using the Boltzmann equation but is very expensive computa- 
tionally within NEGF formalism |51I , should not play an important role in devices depicted 
in Fig. [T] because the mean-free path in graphene is ^ ~ 677 nm at room temperature 1861 . 
That is, the active region of the single-molecule nanojunctions in Fig.[T]or the width of their 
GNR electrodes is much smaller than i at all temperatures below the ambient one which are 
the focus of our study. 



4 Electronic transmission and thermopower in single-molecule nanojunctions with 
GNR electrodes 

The recent fabrication of GNRs with ultrasmooth edges p87"88"891 has opened new avenues 
for highly controllable molecular junctions with a well-defined molecule-electrode contact 
characterized by high transparency, strong directionality and reproducibility. This is due to 
the fact that strong moIecule-GNR %-n coupling makes possible formation of a continuous 
;r-bonded network across GNR and orbitals of ;r-conjugated organic molecules 1901 . Al- 
though GNRs with zigzag 1881 or chiral edges II89I are insulating at very low temperatures 
due to one-dimensional spin-polarized edge states coupled across the width of the nanorib- 
bon, such unusual magnetic ordering and the corresponding band gap is easily destroyed 
above > 10 K I91II92II so that they can be considered as candidates for metallic electrodes. 
In fact, the experimental pursuit of graphene-based single-molecule nanojunctions, where a 
;r-conjugated organic molecule is inserted into the nanogap formed by feedback controlled 
electroburning of few-layer graphene, has commenced very recently [93]. 

Unlike the metallic carbon nanotubes (CNTs), which were employed in earlier exper- 
iments 1941 as electrodes of single-molecule nanojunctions in order to generate ;r-bonded 
network I90L GNRs have planar structure appropriate for aligning and patterning. The early 
experiments 1941 on CNT| molecule |CNT heterojunctions have measured surprisingly small 
conductances for a variety of sandwiched molecules. The first-principles analysis of differ- 
ent setups reveals that this is due to significant twisting forces when molecule is connected 
to CNT via 6-membered rings I90II . Therefore, to keep nearly parallel and in-plane config- 
uration (hydrogen atoms of 18-annulene slightly deviate from the molecular plane) of our 
ZGNR| 18-annulene|ZGNR junction, we use a 5-membered ring I90II to chemically bond 
ZGNR to 18-annulene or to CIO molecule, as illustrated in Fig.[T] 

The high contact transparency makes it possible for evanescent wavefunctions from 
the two ZGNR electrodes to tunnel into the molecular region and meet in the middle of 
it (when the molecule is short enough 1901 ). This is a counterpart of the well-known metal 
induced gap states in metal-semiconductor Schottky junctions. Such effect can induce a 
large peak (i.e., a resonance) in the electronic transmission function near Ep, despite the 
HOMO-LUMO energy gap of the isolated molecule. This phenomenology is confirmed by 
Figs.|2a) and [2c) showing the zero-bias electronic transmission ■%\{E) where the peak 
around E — Ep = is conspicuous for both single-molecule nanojunctions illustrated in 
Fig-El The peak is located far away from HOMO and LUMO levels of the molecules which 
are marked in Figs. 13 a) and (c). 

The carbon atoms of a ring-shaped 18-annulene molecule can be connected to ZGNR 
electrodes in configurations whose Feynman paths for electrons traveling around the ring 
generate either constructive or destructive quantum interference effects imprinted on the 
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E - (eV) E-Ep (eV) 

Fig. 2 (a) Zero-bias electronic transmission .%\{E) for ZGNR| 18-annulene|ZGNR junction; (b) ther- 
mopower at two different temperatures for ZGNR| 18-annulene|ZGNR junction; (c) zero-bias electronic 
transmission for ZGNR|C10|ZGNR junction; and (d) thermopower at two different temperatures for 
ZGNR|C10|ZGNR junction. The vertical dashed lines in panels (a) and (c) mark the position of HOMO 
and LUMO levels of 18-annulene or CIO molecules after they are attached to the electrodes, as extracted 
from the projected density of states for carbon atoms in the center of the junction. 

conductance I46II47I . For example, a ;r-electron at Ef entering the molecule in setup shown 
in Fig.[T]has wavelength kp /2d (d is the spacing between carbon atoms within the molecule), 
so that for the two simplest Feynman paths of length 9d (upper half of the ring) and 9d 
(lower half of the ring) the phase difference is and constructive interference occurs. Note 
that the destructive quantum interference I46II47I would form an additional dip (45] (i.e., 
an antiresonance) within the main transmission peak around E — Ef = in Fig. ^a). The 
effect of such antiresonance on the thermopower S for gold| 18-annulene|gold junctions has 
been analyzed in Ref. (17) as a possible sensitive tool to confirm the effects of quantum 
coherence, even at room temperature, on transport through single-molecule nanojunctions. 

Additionally, the suppression of the hole-like transmission, ^e\{E) — s- for — 1.0 eV < 
E — Ef < 0, avoids unfavorable compensation 1311 of hole-like and electron-like contribu- 
tions to the thermopower. This is due to the symmetry of transverse propagating eigenstates 
in the semi-infinite ZGNR electrodes which we elaborate in more detail in Sec. [5] These 
features in the electronic transmission function yield the thermopower S in Figs.[2tb) and 
|2d) whose maximum value, which is slightly away from E — Ef = at room temperature, 
is an order of magnitude larger than the one measured in large-area graphene 1951 or in 
organic molecules sandwiched by gold electrodes jQj- 
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Comparing Figs.|2ja) and |2b) for ZGNR| 18-annulene|ZGNR junction with Figs.|2c) 
and |2td) for ZGNR|C10|ZGNR junction shows great similarity of the transmission func- 
tion around the Fermi energy E — Ef =0. This is due to the transport through overlapping 
evanescent modes, so that the same features in the electronic transmission function and 
thermopower will persist for any sufficiently short 19011441 organic molecule as long as its 
HOMO and LUMO levels [whose position is marked in Figs. |2ja) and (c)] are far away 
from the Fermi energy of the electrodes. Since these features in the transmission function 
are governed by evanescent and propagating wavefunctions originating in metallic ZGNR 
electrodes, they are impervious to the usual poor estimates of energy gaps and molecu- 
lar level position in DFT and, therefore, do not require more accurate but computationally 
much more expensive NEGF-GW approach 1751 . 



5 Suppression of hole-like transmission by the symmetry of transverse propagating 
states in ZGNR electrodes 

The evanescent mode-induced transmission resonance can also be generated in setups where 
CNT electrodes are attached to a short organic molecul^ via transparent contacts 1901 
l44l However, besides the effect of the same type in our ZGNR| 18-annulene|ZGNR and 
ZGNR|C10|ZGNR junctions, their power factor S^G is further enhanced by completely dif- 
ferent transport properties of the ZGNR conduction band (CB) Bloch eigenstates (exhibiting 
high transmission) and the valence band (VB) Bloch eigenstates (exhibiting low transmis- 
sion). The origin of this asymmetry can be explained by analyzing the spatial symmetry of 
the CB and VB propagating eigenstates. Figure [3t a) shows the DFT-calculated band struc- 
ture of infinite homogeneous ZGNR. The right part of the figure shows isosurface plots of 
the Bloch eigenstates (imaginary part only) in the CB [Fig. [B^b)] and VB [Fig. [S^c)] evalu- 
ated a!Lk = n/2a {a\& the lattice constant of graphene). 

The CB and VB states have different symmetry properties — while the CB state is sym- 
metric with respect to mirroring in the z-axis, the VB state is antisymmetric. Recall that 
molecules in Fig. [T] are coupled to the electrodes via a pentagon in the middle of the trans- 
verse edge of the nanoribbon, and thus form a two-atom linear carbon chain between the 
ZGNR and the rest of the molecule. The state in this linear chain will certainly be symmet- 
ric so that they do not couple with the antisymmetric states in VB, while a large coupling to 
the symmetric states in CB is expected. 

Although the above explanation of the shape of transmission function ■%\{E) around 
E — Ep = Q m terms of the spatial symmetry properties of CB and VB Bloch eigenstates 
is rather intuitive, it may be instructive to show in more detail the connection between the 
eigenstates of the electrodes and the NEGF formalism which often makes no explicit ref- 
erence to such eigenstates |96j|. For this purpose, we consider a simple model depicted in 
Fig. IHa) which shares many features with single-molecule nanojunctions in Fig. [T] The 
model consists of two semi-infinite double tight-binding chains acting as leads that are at- 
tached to a single tight-binding chain in the center of the device. We assume a single ^-orbital 
on each site with on-site energy eo = eV. The hopping parameters are jS = 2 eV for inter- 
chain hoppings and a = 1 eV for all other hoppings. 



^ For an example of the peak in Ta{E) induced by the overlap of evanescent wavefunctions originating 
from two CNT electrodes sandwiching 18-annulene molecule see Supplemental Material of Ref. 1901 . 
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Fig. 3 (a) The band stmcture of 8-ZGNR. Imaginary part of (b) conduction band and (c) valence band Bloch 
eigenstates evaluated aX k = K/2a, as indicated with circles in panel (a). The conduction band state is sym- 
metric with respect to mitToring in the ^-axis, while the valence band state is antisymmetric. 



Thus, the Hamiltonian of the lead unit cell marked by rectangles in Fig.|4la) is given by 



and the coupling matrix between two such unit cells is 

a 



(17) 



Figure |4jb) plots the transmission function of the junction (solid red) and of the infinite 
electrodes (dashed black). Similar to plotted in Fig.|3a) and lUc), we find a com- 

plete blocking of the transmission function for VB states, while the CB states are highly 
transmitting. This asymmetry between quantum transport through VB and CB transverse 
propagating states may again be explained by their symmetry properties. 

In analogy with ZGNR, the VB Bloch eigenstate is antisymmetric while the CB Bloch 
eigenstate is symmetric. The Bloch eigenstates can be written as 



WVB.CB[k) = (pVB.CBS 

V2 \a 



— ikz 



<hB 



(18) 

(19) 



where a = — 1 for VB and a 
expressed as 1971 



+ 1 for CB. We now recall that the lead self-energy can be 



(20) 



using the transfer matrix T. For our purpose it is sufficient to consider only the imaginary 
part of the self-energy because the asymmetry is determined by Fi = — 2Im [X'/^]. The trans- 
fer matrix can be calculated from the complex band structure (721 including both propagat- 
ing and evanescent states. The imaginary part of the transfer matrix is determined solely by 
the transverse propagating states 1721 : 



lm[T{E)]=Y,i»nlv^{e''")i>l 



(21) 
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Left Central Right 

Lead unit cell 




Energy (eV) 

Fig. 4 (a) Simple model of a two-terminal junction consisting of two semi-infinite double tight-binding 
chains, which act as the left and the right lead attached to a single tight-binding chain as the central region. 
The transmission function for this device is shown in panel (b) as the solid red line, while the transmission 
function for the infinite double-chain lead is shown in dashed black for comparison. 



where the sum runs only over propagating states whose transverse part (j)„ is purely real. The 
surface Green function 1961 in the lead-unit cell next to the central single chain is related to 
the self-energy by gL = V^'X'/,V^' . It now follows that the surface Green function has the 
structure 

The self-energy on the single chain due to the left lead is obtained as = V^^gLVic. 
where Vlc is the coupling matrix between the left lead and the central single-chain (there 
is a similar contribution from the right lead as discussed in Sec.[2]i. Since only the first site 
in the chain is connected to the left lead, has only one non-zero element {Si)ii ^ 
in the upper left corner. From Eq. i22i we get that Im[(X'£)ii] oc (1 +a), such that Fi = 
— 2Im [Sl] = for energies in VB, while is non-zero for energies in CB. Therefore, 

since the transmission function of the junction is given by Eq. ([7]i, it follows that it must be 
zero for energies in the VB. 

The analysis presented in this Section demonstrates that the vanishing transmission 
function [Fig. |Hb)] at energies within the VB for the simple model junction, as well as 
vanishing hole-like transmission function [Figs. |2ja) and (c)] in the energy range below 
E-Ef =0 in realistic ZGNR| 18-Annulene|ZGNR or ZGNR|C10|ZGNR junctions, can 
be explained by the symmetry properties of the transverse propagating states in the semi- 
infinite electrodes. This demonstrate the generality of our ZGNR| molecule |ZGNR device 
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Temperature (K) Temperature (K) 

Fig. S (a) The phonon transmission function Sph{oi) and (b) the corresponding thermal conductance Kph for 
an infinite 8-ZGNR and ZGNR| 18-annulene|ZGNR single-molecule nanojunction shown in Fig.^a). (c) The 
phonon transmission function and (d) the coiTesponding thermal conductance ffpi, for an infinite 8-ZGNR and 
ZGNR|C10|ZGNR single-molecule nanojunction shown in Fig.^b). Panels (a) and (b) also compare results 
obtained using NEGF coupled to DPT (via GPAW code (81]) vs. NEGF coupled to Brenner EIP (via GULP 
code (H). 

concept and explains why different molecules can have very similar transmission functions, 
as exemplified by Figs. [2a) and |2c). On the other hand, it also illustrates a weakness in 
the setup because the central molecule must be coupled to one of the ZGNRs exactly in 
the middle if its transverse edge. Nevertheless, these findings open up further possibilities 
to search for other single-molecule nanojunctions which can exploit the symmetry of the 
Bloch states in their electrodes in order to optimize the thermopower 

6 Phonon transmission and thermal conductance in single-molecule nanojunctions 
with GNR electrodes 

Figure m shows the phonon transmission function 5^h(ft>) and the comsponding phonon 
thermal conductance (Cpj, for ZGNR| 18-annulene|ZGNR and ZGNR|C10|ZGNR single-molecule 
nanojunctions computed via the first-principles formalism delineated in Sec. [3] Note that we 
use the term "phonon" here freely to refer to any quantized vibrations in the active region 
of the junction. To understand how the mismatch in vibrational properties of semi-infinite 
ZGNR electrodes and the molecule impedes phonon transport across the junction, we also 
plot ,f5J,h(ft)) [in Figs.Oa) and Oc)] and (fph [in Figs.[5tb) and |5td)] for an infinite homo- 
geneous 8-ZGNR. 

The phonon transmission function for infinite homogeneous 8-ZGNR consists of quan- 
tized steps, as expected for purely ballistic transport of phonons described in harmonic ap- 
proximation |80|. The suppression of phonon transmission by the presence of the molecule 
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generates three times smaller (Cph at room temperature when compared to the thermal con- 
ductance of an infinite 8-ZGNR, as shown in Figs.|5tb) and [Hd). 



6.1 Comparison of DFT and semi-empirical methods for the computation of phonon 
thermal conductance 

The force constant matrix K, discussed in Sec. [3] as an input for NEGF calculations, can 
also be obtained using empirical interatomic potentials. These is widely used methodology 
to treat elastic phonon quantum transport in clear or disordered nanowires l64l|99lfT00ll since 
it is computationally much less expensive that our first-principles approach discussed in 
Sections [3] and [6] One of the standard choices for hydrocarbon systems is the so-called 
Brenner EIP IIIOII . which is often applied to study lattice dynamics and phonon thermal 
transport in carbon nanotubes and graphene 11021 . 

The Brenner EIPs are short range, so they cannot accurately fit the graphene dispersion 
over the entire Brillouin zone (BZ). However, the thermal transport depends more sensitively 
on the accuracy of acoustic phonon frequencies near the zone center where the longitudinal- 
and transverse-acoustic (LA and TA) velocities and the quadratic curvature of the out-of- 
plane acoustic (ZA) branch are determined. Conversely, only weak thermal excitation of the 
optical phonons and acoustic phonons near the BZ boundary occurs around room tempera- 
ture because of the large Debye temperature (~ 2000 K) of graphene. 

The basic steps of NEGF-Brenner-EIP methodology are: initially relax the device geom- 
etry => the force constant between atom 7 in direction a = x,y,x and atom J in direction jS is 
calculated using analytical derivatives, Kj^jp = dU / {dRjadRjp), where U is the total en- 
ergy compute retarded GF in Eq. ( il4l l and phonon transmission function using Eq. ( 1131 . 
Here we employ Brenner EIP as implemented in the GULP code 198111031 . 

To contrast this methodology with the full NEGF-DFT treatment of the phonon ther- 
mal transport in our single-molecule nanojunctions, we add NEGF-Brenner-EIP results for 
^ph{co) and (Cph for the case of ZGNR| 18-annulene|ZGNR junction and infinite homoge- 
neous 8-ZGNR in Figs. [S^a) and [5tb). Although NEGF-Brenner-EIP calculations differs 
very little from NEGF-DFT in the case of infinite homogeneous ZGNR, the NEGF-Brenner- 
EIP methodology yields significantly lower ^hl*!^) and (Cph for ZGNR| 18-annulene|ZGNR 
junction. Thus, using the Brenner-EIP to describe phonon transport would lead to an over- 
estimated ZT. 



7 Thermoelectric figure of merit of single-molecule nanojunctions with GNR 
electrodes 

Combining the first-principles results in Sec.|4]for electron transport and in Sec.[6]for phonon 
transport allows us to obtain the thermoelectric figure of merit ZT. The results shown in 
Fig. |6] demonstrate ZT ~ 0.5 for both junctions at room temperature. For comparison, we 
note that the recent proposal I31II for the single-molecule thermoelectric devices with so- 
phisticated combination of local chemical tuning of the molecular states and usage of semi- 
conducting electrodes has predicted much smaller ZT ~ 0. 1 at room temperature. 

In addition, 0.5 < ZT < 2.5 at E — Ef = —0.02 eV (which can be set by the backgate 
electrode covering the whole device I95II ) in Fig.|6]within the temperature range T = 30-77 K 
is much larger than the value achieved in conventional low-temperature bulk thermoelec- 
tric materials (2l. It is also larger than ZT :^ 1 found in recent studies of low-temperatures 
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Fig. 6 The thermoelectric figure of merit ZT for the two single-molecule nanojunctions shown in Fig. [T] 
(a) ZT vs. energy at two different temperatures for ZGNR| 18-annulene|ZGNR junction; (b) ZT vs. tem- 
perature at three different energies for ZGNR| 18-annulene|ZGNR junction; (c) ZT vs. energy at two differ- 
ent temperatures for ZGNR|C10|ZGNR junction; and (d) ZT vs. temperature at three different energies for 
ZGNR|C10|ZGNR junction. 



(T < 100 K) nanoscale thermoelectrics, such as the Kondo insulator nanowires f23|. Thus, 
our single-molecule nanojunctions could be suitable for high-performance thermoelectric 
cooling at low temperatures. 



8 Summary and future prospects 



We provided an overview of NEGF-DFT methodology to treat independent electron and 
phonon transport in nanoscale thermoelectric devices, where we focused on the nascent 
subfield of single-molecule nanojunctions. The examples of nanojunctions we considered 
consist of a ;r-conjugated organic molecule attached to metallic ZGNR electrodes via strong 
covalent bond which creates continuous ;r-bonded network across the device. Although de- 
vices in Fig.[T]look futuristic at first sight, presently available nanofabrication technologies 
have already demonstrated feasibility of single-molecule nanojunctions with graphene elec- 
trodes II93I . The highly transparent GNR| molecule contact allows for evanescent wavefunc- 
tions to penetrate from the ZGNR electrodes into the HOMO-LUMO gap of sufficiently 
short molecule, thereby generating a transmission resonance as a favorable ingredient to op- 
timize the electronic thermopower. While such resonance could be achieved by using other 
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types of carbon-based metallic electrodes such as metallic CNTs, the usage of ZGNRs also 
brings peculiar spatial symmetry properties of transverse propagating states which block 
transmission of electrons in a range of energies below the Fermi energy of the device. This 
increases the thermopower where, otherwise, contributions from hole-like and electron-like 
transmission appear with different signs in Eq. dlOb and thus partially cancel each other. 
Finally, the mismatch in vibrational properties of the semi-infinite ZGNR electrode and the 
molecule acts to reduce the phonon thermal conductance across the junction by a factor 
of about three when compared to infinite homogeneous ZGNRs. The combination of these 
three features leads to ZT ~ 0.5 at room temperature and 0.5 < TT < 2.5 below liquid nitro- 
gen temperature for the junctions in Fig.[T] These values of ZT turn out to be much higher 
than those found in other recent first-principles studies of single-molecule thermoelectric 
devices j28]|3T] |33l. 

While the thermoelectric figure of merit ZT we obtained at room temperature is still 
low, we anticipate that much higher ZT could be achieved by testing different types of 
molecules to reduce (Cpj, further. Essentially, the single-molecule nanojunctions proposed 
in Fig.[T]should be viewed as representatives of a new class of nanoscale thermoelectric de- 
vices where the power factor S^G is already optimized by the usage of GNR electrodes that 
generate molecular-level-independent transmission resonance and where the spatial sym- 
metry of transverse propagating eigenstates in GNR electrodes lifts the compensation of 
hole-like and electron-like contributions to S. 



8.1 Future prospects 

We also offer a brief survey of possible future directions and related challenges. Our dis- 
cussion of quantum electronic and phononic transport in realistic single-molecule nano- 
junctions via NEGF-DFT formalism has been confined to non-interacting electrons and 
phonons which propagate independently of each other. Since single-molecule nanojunc- 
tions depicted in Fig. [T] have highly transparent contacts due to strong molecule-electrode 
coupling, we have not considered Coulomb blockade effects Il34ll36ll35l that would become 
relevant in molecules weakly coupled to the electrodes. Their understanding (such as simul- 
taneous treatment of Coulomb blockade and coherent tunneling transport regimes II34I ) to 
correctly capture thermoelectric properties in the presence of many-body interactions is in 
its infancy. 

Although the electron-phonon interaction is the dominant mechanism of inelastic scat- 
tering in single-molecule nanojunctions |55 104 49,50|, whose effect on electronic cur- 
rent I55lll04ll50lll05l and its noise |59] has attracted considerable attention, the study of 
electron and phonon quantum transport coupled by such interaction is extremely rare. Also, 
while the techniques have been developed to take into account anharmonicity into NEGF 
calculations of phonon transport through single-molecule nanojunctions fSll, the treatment 
of phonon-phonon scattering via such fully quantum-mechanical approach combined with 
first-principles methods to capture atomistic structure of realistic junctions is lacking. 

For example, there is only a handful of recent papers I33II48II 106)1 107 IfTOS " 1091 where 
this problem of coupled electron and phonon transport has been tackled using NEGF applied 
to simple phenomenological models 14811 10611 10711 10811 1091 , or NEGF-DFT with simplifica- 
tions in considering interacting electrons and noninteracting phonons ||33I. Nevertheless, 
these attempts hint 133 ,481 at possibly strong effect of electron-phonon interaction on the 
thermopower 5, as well as that ZT of single-molecule nanojunctions can be enhanced by 
them. In fact, electron-phonon interaction is largely an unexplored parameter in the quest 
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for efficient thermoelectric materials, as highlighted by the recent experiments on superlat- 
tices lUOI . 

The main technical issue in evaluating even the lowest order Feynman diagrams (Hartree 
and Fock for electrons 149115011 and polarization bubble for phonons 11091 ) for the retarded 
and lesser nonequilibrium self-energies 1581 is integration in energy which has to be done 
repeatedly in order to achieve self-consistency that ensures conservation of charge and heat 
currents II106I . The successful solution to this extremely computationally intensive 1491 
problem will open a pathway to understand phonon drag (interchange of momentum be- 
tween acoustic phonons and electrons that generates additional contribution to the ther- 
mopower which has thus far been studied only via coupled electron-phonon Boltzmann 
equations Hill ) and electron drag (where phonons are dragged by electrons from low T 
region into high T region 11091 ) effects in realistic single-molecule thermoelectric devices 
using first-principles quantum transport simulations. 
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